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Abstract 

Starting from Wilson's action, we calculate strong coupling series for the 
Polyakov loop susceptibility in lattice gauge theories for various small N T in 
the thermodynamic limit. Analysing the series with Pade approximants, we es- 
timate critical couplings and exponents for the deconfinement phase transition. 
For SU(2) pure gauge theory our results agree with those from Monte-Carlo sim- 
ulations within errors, which for the coarser N T = 1,2 lattices are at the percent 
level. For QCD we include dynamical fermions via a hopping parameter expan- 
sion. On a N r = 1 lattice with Nf = 1,2,3, we locate the second order critical 
point where the deconfinement transition turns into a crossover. We furthermore 
determine the behaviour of the critical parameters with finite chemical potential 
and find the first order region to shrink with growing fx. Our series moreover 
correctly reflects the known Z(N) transition at imaginary chemical potential. 



1 Introduction 



Lattice Monte Carlo studies of the QCD phase diagram at finite temperature are still 
difficult for realistic quark masses, and at finite baryon density are beset by the sign- 
problem. Despite the progress made in the last few years in circumventing these obsta- 
cles, such calculations still suffer severe restrictions. In particular, presently available 
methods are reliable only for small quark chemical potentials, /x < T pQ . These diffi- 
culties motivate the search for alternative ways to learn about the phase diagram. A 
popular choice is to consider lattice QCD in the strong coupling limit and to study 
its phase diagram either by analytic mean field methods [2] or indeed by Monte Carlo 
evaluation of the strongly coupled theory [3], [I], since the sign problem in this case is 
much milder. In the early days of lattice gauge theory, strong coupling expansions for 
zero temperature Yang-Mills theory have led to some analytical insights into field the- 
ories on a lattice ([5], [6] and [7] for a review). They have also been used to investigate 
finite temperature effects [E]-[I2] with reasonable qualitative predictions, but mostly 
using the crude approximation of neglecting spatial plaquettes altogether, cf. [H] and 
references therein for a review of early investigations. However, the strong coupling 
limit is far from the physical continuum theory and the lessons learned in this way are 
qualitative at best. There have also been attempts to go beyond the strong coupling 
limit [15]- [19], mostly in mean field theory. 

In this paper we address the question whether it is possible to make predictions for 
the location and nature of the deconfinement phase transition more quantitative by 
taking corrections into account, and we can answer in the affirmative. In previous 
work [20J we have already shown for SU(2) Yang-Mills theory that it is quite possible 
to include spatial plaquettes and to obtain series of several orders. For lattices with 
temporal extent N T = 1 — 4, this lead to satisfactory quantitative results for the 
equation of state up to the phase transition region. However, while the determination 
of the critical parameters was consistent with Monte Carlo results and universality 
arguments, it remained rather imprecise quantitatively. Here we considerably improve 
on this by using the Polyakov loop susceptibility as an observable, rather than the free 
energy. At a second order phase transition in infinite volume, this observable develops 
a singularity which is well modelled by Pade approximants to its series expansion, thus 
allowing to extract the critical coupling and exponent in satisfactory agreement with 
numerical results from |21j,[22j,[23j. 

After a successful test of our method for SU{2) Yang-Mills theory, we study the 
case of SU(3) QCD with Wilson quarks by a combined strong coupling and hopping 
parameter expansion, which converges for sufficiently heavy quarks. In the infinite 
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Figure 1: Schematic overview of the Nf = 2 + 1 QCD phase diagram. Here we compute 
the critical line for heavy quarks. 




quark mass limit, QCD is known to display a first order phase transition which weakens 
as quark masses are lowered, see Fig. [IJ The location of the critical quark mass, where 
the transition disappears at a second order critical point, has been studied for Nf = 1 
by a combination of numerical simulations of QCD and the 3d 3-state Potts model 
[2U [22]> which is the appropriate effective theory for the critical transition. These 
works also established the 3d Ising universality of the boundary line. Our methods 
allow for a determination of the critical quark mass for Nf = 1,2,3, as well as the 
dependence of the critical quark mass on quark chemical potential. We find that the 
critical mass grows with increasing chemical potential, in accord with a Monte Carlo 
study in the Potts model with finite chemical potential [25]. However, our result is the 
first based on full QCD beyond mean field theory. 



2 Phase transitions from strong coupling series 
2.1 Notation and formalism 

We work on a (3 + l)-dimensional hypercubic lattice with lattice spacing a and infinite 
spatial volume. The temporal lattice extent N T is kept finite, which in combination 
with (anti-)periodic boundary conditions for (fermionic) bosonic fields generates a non- 
vanishing physical temperature, 
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Using Wilson's gauge action for SU(N), the partition function reads 
Z = J [dU}exp(-S g ) = J [dU] JJexp ^ (tr E/p + tr 



(2) 



In order to locate the phase transition we consider the Polyakov loop susceptibility 



XL = V(V)-(L> 2 ). 



where we have defined the Polyakov loop L x and its spatial average L as 



L x = trW x = tr Y[U (x,t), 

T=l 

L = ^E( L - + L x) SU(N>3) 

X 

L = ^$>x SU{2). 

X 

If we couple the Polyakov loop to an external source J in the action^, 

- 5 ( J ) = £f E ( tr u p + tr + J E ( L * + L i 

p X 

we can express the susceptibility as 

1 d 2 

XL= vdT 2lliZ{J) 



(3) 



(4) 



(5) 



(6) 



(7) 



j=o 



In order to obtain a strong coupling series for Eq. (j7J), we expand the partition 
function 



Z(J) = j [dU] J] exp (A (tr t/ p + tr £/ p t) ) J] exp ( J (trW x + trW x 



(8) 



in terms of characters 



2iV 



exp | — I trZ/p + trZ/jJ 



1 + ^<i r a r Xr(U p 



exp f J (trW x + trW x 



c (J) 



l + ^6 r (J) Xr (M/ x 



(9) 



(10) 



1 For SU (3) we have chosen L as the real part of the Polyakov loop, so we get a real action with 
only one real source J 



Figure 2: Examples for N T = 4. Left: Two Polyakov loops and N T plaquettes wind 
around the temporal dimension. Right: The first correction with additional plaquettes. 



We have neglected the prefactor c (/3) of the trivial representation because it does not 
depend on J and vanishes in Eq. (J7J). As usual in strong coupling expansions, we will 
use the coefficient u = aj of the fundamental representation as expansion parameter. 

In the case of 577(2) we have no complex conjugate representations and the partition 
function reads 



Z{J) 



[dU] 



JJexp f^trC/p) 



j Jcxp ( JtrW : 



(11) 



Furthermore, there are closed form expressions for the expansion coefficients in this 
case, 



a, (0) 

W) 
c (J) 



09) 



d 3 h 3+ i(2J) 
h(2J) 

h(2J) 
J ' 



(12) 
(13) 
(14) 



Applying a cluster expansion as described in [27], the logarithm of the partition 
function can be represented as a sum of graphs $(C) 



LnZ(J)= In c (J)+^$(C). 



(15) 



The sum is over all clusters C of connected polymers, see [20j for details. The leading 
order, i.e. the strong coupling limit (3 = 0, is obtained by neglecting all graphs giving 
the trivial result 
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SU(N > 3). 



(16) 



In case of SU (N > 3) the factor of 2 accounts for both fundamental representations. 



The first graph with a non-trivial ^-dependence is shown in Fig. [2J together with the 
leading correction graph. The left and the right boundaries of these graphs are meant 
to be identified due to the periodic boundary conditions. To calculate the contribution 
of these graphs, we employ the group integration formula 

J dU X r(U) Xs m = ^, (17) 

and get 

$ = 3u N " (b 1/2 (J)) 2 SU{2), 

$ = 6 u Nt {b f {J)f SU{N > 3), (18) 
for the leading non-trivial order and 

$! = \2N T u N ^ A (b 1/2 {J)f SU(2), 

$i = 2AN T u NT+A (b f (J)f SU{N>3). (19) 
for the first correction. 

2.2 Graphical expansion 

From Eq. (j7]) it is obvious that we solely have to take into account graphs which 
contribute to order J 2 . This means that only graphs with two Polyakov loops in 
the fundamental or one loop in the adjoint representation are allowed. For the first 
possibility, the loops have to be on different lattice sites. The generation of contributing 
graphs is not uniform and we distinguish between small, intermediate and large N T . 

Large N T : Large N T receive only contributions from nearest-neighbour Polyakov 
loops as shown in Fig. [2J and corrections from adding plaquettes. Of course, this 
statement is only true for large enough N T if we calculate to some fixed order in u. 

Small iV T : The smallest possible N T is 1. Typical graphs are shown in Fig. [3j These 
graphs are meant to be spatial projections of graphs like in Fig. [2] (left). In higher 
orders we get contributions from additional spatial plaquettes, e.g. by filling the cross- 
section of the self-avoiding polygons, but these contributions are small compared to 
the increasing number of self-avoiding walks. 

Intermediate N T : For intermediate N T (=2,3,4...) we have to take into account 
graphs of both types. There are also some other corrections as shown in Fig. HI Thus 
these N T are the most labour-intensive ones. 
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Figure 3: Left: Self avoiding walks with two fundamental Polyakov loops. Right: Self 
avoiding polygons with one adjoint or two fundamental Polyakov loop. 




Figure 4: Examples of corrections to self avoiding walks of length L = 2 and N T = 4. 



2.3 Series analysis and phase transitions 

Strong coupling expansions have a finite radius of convergence. Expanding about (3 — 0, 
a true deconfinement phase transition at some critical value of /3 C clearly represents an 
upper bound on the convergence radius, i.e. strong coupling analyses are limited to 
the confined region. Nevertheless, knowledge of the series coefficients allows us to 
estimate the location of a singularity along the real /3-axis. Our analysis is well suited 
to detect second order phase transitions and exploits the fact that the Polyakov loop 
susceptibility in this case diverges with a critical exponent. Near a critical coupling 
the Polyakov loop susceptibility and its logarithmic derivative behave like 



XL 



zy„)^i„( Xl) 
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(20) 



(u c — up ' 1 du (u c — u) ' 

From our series expansions we know D x (u) as a polynomial in u and can model its 
pole-like singularity by Pade approximants 

t L > M \\ u ) = i ,»..., 1 I 77m ' ( 21 ) 



l + b 1 u + --- + b M u M ' 

In order to uniquely determine the coefficients a*, bi, it is necessary to have L + M < N, 
if iV represents the highest available order of the expansion. In this way a [L, M] 
approximant is correct up to but not including 0{u L+M+1 ) and larger approximants 
represent more expansion coefficients than smaller ones. In particular, the critical 
coupling u c is given as the real positive zero of the denominator closest to the origin, 
the critical exponent 7 is obtained from the corresponding residuum. 

The analysis can be made more powerful if either independent results for the critical 
couplings are available, or the universality class of the transition is known. In the first 
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case, it is possible to get better estimates for the critical exponents via 

(u c -u)D x (u) = 7 + 0(m c -m), (22) 

which is the more precise the better we know u c . In order to do so, we calculate Pade 
approximants to the series expansion of (u c — u)D x (u) and evaluate them at the known 
critical coupling u c . In the same way we can use a known value 7 to get more accurate 
estimates for the critical coupling 

(Xl)- ~ • (23) 

[U c - U) 

Here we calculate Pade approximants to (xl)~ and solve for zeros of the denominator, 
as this quantity has a simple pole at the critical coupling. For a more detailed discussion 
of these topics, see [28J. 



3 SU{2) Yang-Mills 
3.1 Results for the series 

We first apply our analysis method to SU{2) pure gauge theory, where we have rea- 
sonably long series for N T = 1 — 4 and where accurate Monte Carlo data are available 
for comparison. We obtain the following strong coupling series for xl{N t ,u): 

52982 

Xl(1,m) = 1 + 6 m + 30 u 2 + 150m 3 + 738u 4 + 3622 u 5 + m 6 + 

773434 7 11239612 „ ^ , qN 
+ —^~ u + 27 + ° ( u ) ' 

30,(2, u) = 1 + 6 m 2 + 30 m 4 + 222 m 6 + 1218 m 8 + ^^L u 10 + O (u 12 ) , 

Xl(3,m) = 1 + 6m 3 + 30m 6 + 72m 7 + 72m 8 + 78m 9 + 576m 10 + 1776m 11 + 
+ 1770 u 12 + (u 13 ) , 

Xl(4,m) = 1 + 6n 4 + 126n 8 + 48 m 10 + 2830 m 12 + ^^n 14 + 0(n 16 ). (24) 

135 

N T = 1 corresponds to the largest lattice spacing and thus to the largest bare coup- 
ling at the deconfinement transition. Hence, our series shows the best convergence 
behaviour in this case. 
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Pade 


u c 




7 


[5,2] 


0.21055 


0.86825 


1.167 


[3,3] 
[4,2] 
[2,4] 


0.20967 
0.20957 
0.20987 


0.86439 
0.86396 
0.86527 


1.138 
1.146 
1.135 


[3,2] 


0.20927 


0.86264 


1.126 


[2,2] 


0.20820 


0.85796 


1.102 


Mean 


0.2095(12) 


0.864(5) 


1.14(3) 



Table 1: Critical coupling and exponent for N T = 1, estimated from different Pade 
approximants. 

3.2 The critical parameters 

We consider only Pades [L, M] with L, M > 2 in order to have large enough polynomials 
in both numerator and denominator. So-called defective approximants with an adjacent 
zero-pole pair, indicated by a small residuum (we chose Res < 0.003 to be defective), 
are also ignored. By doing so we obtain estimates for the critical parameters as shown 
in Table [TJ In the last line we have averaged over the estimates from different Pades 
in order to quantify the systematic error associated with the choice of a particular 
approximant. Note that the quoted error is estimated as (/3™ ax — /3™ m )/2, and similarly 
for the exponent. It is only due to the scatter in the singularity structure of different 
Pade approximants and does not include the error from the truncation of the series, 
i.e. it likely underestimates the true error. 

Our result for a N T = 1 lattice then is (3 C = 0.864(5), 7 = 1.14(3). There exist 
different values for the critical coupling from Monte Carlo simulations in the literature: 
f3 c = 0.8730(2) from [21] and the more recent & = 0.85997(10) or 0.86226(6) from [23]. 
At first sight our results appear to favour one of the latter, but with a critical exponent 
deviating about 10% from universality. While the Pade approximants accumulate 
a pole and thus definitely predict a second order phase transition, the residuum is 
somewhat below the value 7/ = 1.2373(2) [2H] for a 3d Ising transition. Note however 
the upward trend in the critical coupling as well as in the exponent with increasing 
order L + M of the approximants in Table [TJ This indicates that the results are not 
yet fully stable and the exponent should reach the Ising value with longer series. 

It is apparent that we have gained considerable accuracy compared to previous 
work using the free energy and its derivatives as observables [20], which gives /3 C = 
0.92(15), a = 0.063(38) (a = 0.1096(5) for 3d Ising [23]). The reason for this improve- 
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Pade 


u c 


0c 


[6,2] 


0.21221 


0.87553 


[4,3] 
[2,5] 


0.21159 
0.21138 


0.87281 
0.87189 


[3,3] 
[2,4] 
[4,2] 


0.21229 
0.21238 
0.21279 


0.87588 
0.87628 
0.87808 


[3,2] 
[2,3] 


0.20986 
0.21464 


0.86523 
0.88621 


[2,2] 


0.21495 


0.88757 



Pade 


7i 


72 


[3,4] 


1.1250 


1.2378 


[4,3] 


1.1244 


1.2331 


[2,5] 


1.1246 


1.2157 


[5,2] 


1.1244 


1.2208 


[3,3] 


1.1225 


1.2579 


[2,4] 


1.1236 


1.2661 


[4,2] 


1.1244 


1.2950 


[3,2] 


1.1240 


1.2308 


[2,3] 


1.1238 


1.2215 



Table 2: Biased critical couplings and exponents for N T — 1. 



ment is twofold: the Polyakov loop susceptibility permits an easier evaluation of more 
coefncents, e.g. by featuring both even and odd powers of u, and the series itself comes 
with only positive coefficients and is better behaved than that for the free energy. 

It is now interesting to explore how one can combine Monte Carlo and series results. 
Thus we consider biased estimates, which should be more accurate. Using 7/ = 1.237, 
we get the results shown in TableEJ We calculate the average of C to be C = 0.877(11) 
using all Pades and (3 C = 0.875(3) using only the three highest orders which behave 
more smoothly. Despite the fact that the total error is underestimated, both estimates 
are consistent with the Monte Carlo result of [21] . 

In order to obtain the biased critical exponent, we used the values 0% = 0.86226 and 
02 — 0.873. The former gives a mean critical exponent of 77 = 1.124(1) and the latter 
7^ = 1.24(4). Although the first result is much more stable, it is the second one which 
is consistent with universality. Hence we conclude that it is the value C = 0.8730(2) 
of [21], which is supported by our series expansions. 

For intermediate N T our results become less precise the larger we choose N T . This is 
to be expected since C grows on finer lattices and we are leaving the strong coupling 
regime. Thus we only give our biased estimates, using the Monte Carlo results C (N T = 
2) = 1.87348, C {N T = 3) = 2.1768, C (N T = 4) = 2.2993 [22], [23] and 7/ = 1.237. We 
summarise our results in Table [3] and observe that the predicted quantities are fully 
consistent with Monte Carlo results and universality. 
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N T 


7/ 


7 


No. of Pades 




Jc 


No. of Pades 


2 


1.237 


1.21(2) 


5 


1.87348(2) 


1.87(1) 


6 


3 


1.237 


1.29(18) 


2 


2.1768(30) 


2.13(4) 


6 


4 


1.237 


1.22(20) 


3 


2.2993(3) 


2.23(11) 


7 



Table 3: Comparison of our findings, 7, /3 C , with the values from universality and 
simulations. 



4 SU{3) and QCD 

For SU (3) Yang-Mills theory, there is a first order phase transition, i.e. the correlation 
length remains finite even at the critical temperature, spoiling our analysis method 
which requires scaling behaviour. We therefore introduce heavy dynamical quarks, 
which have a critical mass m c where the transition turns second order, cf. Fig. [0 It is 
this point which we now try to locate. 

4.1 Combined strong coupling and hopping expansion 

We introduce dynamical quarks in leading order hopping parameter expansion (see [27] 
for details), where the quark part of the action reads 

* = £jf"P]', , = ^, (25) 

I 

m is the quark mass and M[[7] the quark hopping matrix 

M[U] yx = <W(1 + 7m)*V (26) 

The sum in Eq. ( 125]) extends over all closed paths on the lattice. For small temporal 
lattice sizes and finite temperature the leading order hopping expansion term is just 
the Polyakov loop tiW^. Chemical potential is introduced in the usual way as factors 
exp(±/i) to the temporal link variables [30J. The effective quark part of the action for 
small temporal lattice extents then reads 

- S e q ff = [%)e M trW x + h(K)e-^trWl 

X 

where the relative minus sign compared to Eq. (^o§ is due to the antiperiodic boundary 
conditions for fermions. The parameter h(n) depends on the hopping parameter k and 



(27) 
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the number of degenerate quark flavours Nf via 

h(n) = 2N f (2K) NT . (28) 

For 7^ the action is complex. 

Introducing Polyakov loop source terms J and putting everything together, we obtain 
the following partition function 



[dU] exp 



+ { \- hdl + J ] tr w * + [ /ie_M + J ] tr w * } 



(29) 



Now we can proceed in the same way as in case of 577(2) and rewrite the partition 
function with two different character expansions and omit the factor cq(/3) 



[dU]H 



1 + ^2d r a r ([3)xr(U p 



X 



Yl c Q (h, [i, J) 



(30) 



f ^2b r (h,fi, J)Xr(W x 

Note that for non- vanishing chemical potential /i the expansion parameters b r (/i, /i, J) 
are different for complex conjugate representations r and f and related via 

b r (h, /i, J) = b f (/i, — /i, J) . 

The expansion coefficients itself can be expressed as series expansions, e.g. 

— /? + ..., 
18 

1 + {he* + J) (he- 11 + J) + . . . , 
he" + J + . . . . 



(31) 



u = a/(,5) 

Co(/i,//, J) 

b f (h,n,J) 



(32) 



In order to get the proper series expansion we now have to draw all possible diagrams 
to a given order in u and the number of Polyakov loops /. In contrast to our SU(2) 
calculation we have to take into account not only the graphs with two Polyakov loop 
source terms, but all graphs with contributions of the order J 2 h m , since these give 
finite results after differentiating twice with respect to J and setting J = 0. Some 
examples of graphs for the case N T = 1 are given in Fig. |5j An important fact is that 
for a given order in u, there is only a finite number of graphs. For the order u n , we 
can have only graphs fulfilling / < 2n, as inspection of the series shows. Additional 
terms in the hopping expansion will modify the 6 r 's and rapidly increase the number 
of relevant graphs. 
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Figure 5: Examples of terms with a larger number of Polyakov loop source terms. 



4.2 Result for the series 



We have derived the series expansion of the Polyakov loop susceptibility up to orders 
u n h m , with n + m < 6. The result for the N T — 1 series, arranged in increasing orders 
)f u. reads 



of u, reads 

XL(u,h) = 



1 + ch 



--c 3 + -c 



+ <l! c5 + r 3 



j v 3 3 24 y 

l5 /28 6 7 4 7 2 119 » ., 
hr + — c 6 - -c 4 - — c 2 + - 



-C U S | _ C 6 -- C 4 -— . 

8 y V 15 5 120 720 



+ 



6 + lBch + (6c 2 + 3) h 2 + (-40c 3 + 15c) h s 



66c 2 - + (-^c 5 - 8c 3 + 6c) /i 5 



+ ( -90c 4 + 

30 + 180c/i + (144c 2 + 72) h 2 + (-760c 3 + 285c) h 
8985 



u 



5985 4 

+ 



c 2 - 



+ 
+ 
+ 

+ (-6856c 3 + 2571c) h 
+ 



u 



/4113 , 

150 + 1470c/i + I — c 2 



V 2 



4113 \ 



'^o. ,n„m , /H31747 2 1088547 
786 + 10752c/i + ( — — — c 2 + 



V 32 



64 



+ 



4011 + 73521c/i 



5 !52247 6 
u H — u b , 



(33) 



where we used the abbreviation c = cosh(/x). Since the only yU-dependence appears 
in cosh(/i) terms, one can immediately see that the Polyakov loop susceptibility is 
invariant under fj, -H- —/i as it should be according to the charge conjugation symmetry 
of QCD. 



4.3 Critical point for fi = 

In order to locate the critical point (f3 c , k c )(/j,) we have to adjust our analysis methods 
to multiple variables. For a given number of flavours, the schematic phase diagram is 
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0.18 




0.12 ' — 1 1 1 1 1 1 1 1 

U 0.6 0.8 1 1.2 1.4 1.6 1.8 2 



n 

Figure 6: Left: Schematic phase diagram in the variables (u, h) of the series Eq. (I3"3"|) . 
Right: Poles of different Pade approximants (red: [1,3], blue: [0,2], green:[0,3], pur- 
ple: [1:2]) accumulate for a slope parameter n = 0.730 at the critical point u c = 0.126. 



shown in Fig.[6j Our expansion is performed about (u, h) = (0, 0) and we now approach 
the critical point along some straight line starting at the origin, 

u — n ■ t, 

h = --t, (34) 
n 

whose slope is tuned by n and t parametrises the distance from the critical point. We 
then expect a scaling behaviour 

with some critical exponent A. The critical point is in the 3d Ising universality class, 
but we do not know the scaling fields and variables. The Polyakov loop in general 
mixes contributions from both energy-like and magnetic field-like variables. When 
approaching the critical point along a straight line from the origin, the larger exponent 
will dominate and in case of 3d Ising universality this is 7. 

Our strategy now is to vary n and calculate DLog-Pades and singularities for each 
value. If the axis misses the critical point, there is no scaling behaviour and we expect 
any real poles in t, and hence in u, of the Pade approximants to be widely scattered. 
As the critical point is approached, these poles accumulate in a narrow window, as in 
our SU(2) study, cf. Fig. [6j Thus, we estimate t c as the mean value over all Pade 
singularities evaluated at that n for which its standard deviation is minimal. As an 
error estimate we take the larger one of those two values, where the standard deviation 
is 1.5 times its minimum. With this method we find for /1 = 

n = 0.730(16), £ = 0.172(4), A = 1.03(3). (36) 



13 



Let us note that the standard deviations for both t and A reach their minimum for the 
same slope parameter n. The critical exponent obtained with this method again under- 
estimates 7 compared to the Ising exponent. Judging from our experience with 577(2), 
we associate this with the truncation of the strong coupling series. From Eq. ( |34l) we 
obtain 

u c = 0.126(1), h c = 0.236(11). (37) 

To get more accurate results, we now employ biasing with the 3d Ising exponent 
7 = 1.237, which leads to the improved values 

u c = 0.131(1) -> p c = 2.03(2) h c = 0.249(13). (38) 

Note that, to leading order in the hopping expansion, fl c does not depend on Nf. 
Inverting Eq. (128]) in case of \x = we find 

-!fJrf (39) 



N f 


= 1 : 


K c 


N f 


= 2 : 


K c 


N f 


= 3 : 


K c 



2 \2N ff 

This leads to the following results for K c (Nf) on a N T = 1 lattice 

0.062(4), 
0.031(2), 

0.021(1). (40) 

For these small values our leading order in the hopping expansion should be an excellent 
approximation. This further justifies use of the relation [9] 

k = \ e"™, (41) 
which is valid for heavy quarks, to obtain the critical quark masses as 



Nf 


= 1 : 


m c /T 


= 2.08(7), 


Nf 


= 2 : 


m c /T 


= 2.78(7), 


Nf 


= 3 : 


m c /T 


= 3.17(10) 



(42) 

The relative size is consistent with qualitative expectations. Since the presence of finite 
mass quarks weakens the first order transition, more flavours of quarks should have a 
stronger effect and thus a larger critical quark mass. 
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4.4 Critical point for /i 



Next we turn on a chemical potential. Eq. (133]) has already been obtained for finite 
chemical potential, all we have to do is to repeat the same steps to search for singu- 
larities as described in the previous analysis. Of particular interest is the movement of 
the critical masses m c with /i. We have calculated this function for the case of Nf = 3 



Figure 7: The critical mass m c (fx 




1 as a function of real and imag- 



inary chemical potential. Error bars are obtained with the same method as described 
before Eq. ( 1361) . 



at several points, Fig. [7J The critical masses grow with real fi, i.e. the first order region 
is shrinking. For small fi we performed a fit to a low order polynomial to get a rough 
picture of the behaviour of the signs of the different coefficients, 



m c (n 2 ) = 3.18 + 0.94(1)// - 0.34(1)// + 0.037(18)/i 6 + 



(43) 



where the errors are those of the fitting procedure. The shrinking of the first order 
region, and in particular the positive curvature of the critical surface as well as the 
alternating signs are in accord with the findings of a Monte Carlo investigation of 
the Potts model [26J. To leading order hopping expansion different Nf only shift the 
constant term in Eq. fj43|) . 

The alternating signs indicate a convergence limiting singularity on the negative // 
axis, i.e. at imaginary chemical potential. Thus it is interesting to continue m c (/i 2 ) also 
to negative values of the argument, by setting /i — > In Eq. ( 1331) . this means that 
c now abbreviates cos(/a) instead of cosh(/i). The corresponding curve is also shown 
in Fig. [7J We observe a minimum at about fi 2 ~ —0.85 and rapidly increasing errors 
with more negative fi 2 . We interpret this as the point where the singularity is located. 
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This is fully consistent with the Roberge- Weiss Z(N) transition point in the imaginary 
direction, which is known exactly to be at fi 2 = — 7r 2 /9 ~ — 1.1 |31j . 

5 Conclusions 

We have explored the feasibility of calculating the critical couplings and exponents of 
the deconfmement transition of lattice QCD by means of analytically computed strong 
coupling series, using the Polyakov loop susceptibility as an observable. Our treatment 
is purely analytical and goes beyond mean field theory. Tables [21 [31 Eq. ( 14"2"j) and Fig. [7] 
contain our main results. For SU(2) Yang-Mills theory on iV T = 1 — 4 lattices, our 
results are fully consistent with numerical values from Monte Carlo simulations. In 
the cases N T = 1,2 we were able to reproduce these values nearly exactly with an 
accuracy of a few percent, for N T = 3,4 the results agree within increasing error bars. 
For N T > 4 the strong coupling series to the computed length becomes unpredictive, 
at least for our observable. A similar conclusion was drawn in [T5] . 

In case of SU(3) with heavy quarks we performed a strong coupling expansion in 
the effective action to leading order hopping expansion. Analysis of the N T = 1 series 
allowed to extract the critical quark mass m c (/i), for which the first order deconfmement 
transition goes critical before turning into a crossover. We find the first order transition 
region to shrink with increasing chemical potential. This is consistent with the findings 
from mean field theory [11] and a Monte Carlo simulation of the effective theory with 
the same global symmetries, the 3-state Potts model in 3 dimensions [26] . The latter is 
the appropriate effective model also for continuum QCD. Since our series furthermore 
correctly reflects the presence of the Z(N) transition in the direction of imaginary 
chemical potential, this would suggest that the qualitative phase structure is correctly 
represented on lattices as coarse as N T = 1. 

We conclude that strong coupling expansions are able to provide qualitative and, 
on coarse lattices, quantitative information about the phase structure of lattice QCD 
beyond mean field theory and the strong coupling limit, which easily extend to finite 
baryon density. While the location of phase transitions in the parameter space is subject 
to renormalisation in the continuum limit, the presence of critical lines or surfaces is 
guaranteed by universality to survive in the continuum limit. This strongly motivates 
further studies to extend our analyses to finer lattices and to the light quark regime. 
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